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Abstract 


We sequenced the mitochondrial (mt) ND5 gene of 100 specimens of 
Eira barbara (Mustelidae, Carnivora). The samples represented six out of 
the seven putative morphological subspecies recognized for this 
Mustelidae species {E. b. inserta, E. b. sinuensis, E. b. poliocephala, E. b. 
peruana, E. b. madeirensis, and E. b. barbara) throughout Panama, 
Colombia, Venezuela, French Guiana, Brazil, Ecuador, Peru, Bolivia, 
Paraguay, and Argentina. The main results show that the genetic diversity 
levels for the overall samples and within each one of the aforementioned 
putative taxa were very high. The phylogenetic analyses showed that the 
ancestor of the Central and South-American E. barbara originated during 
the Miocene or Pliocene (6.3-4 millions of years ago, MYA). 
Furthermore, the ancestors of some geographical groups, (we detected at 
least four) originated during the Pliocene (3.7-2.5 MYA). These four 
groups (or lineages) were placed in the Cesar-Antioquia Departments 
(northern Colombia), Bolivia and northwestern Argentina, northern- 
central Peru, and in the trans-Andean area of Ecuador. However, during 
the Pleistocene, this species experienced a strong population expansion 
and many haplotypes expanded their geographical distributions. They 
became superimposed on the geographical areas of older geographical 
groups that originally differentiated during the Pliocene. Until new 
molecular studies are completed, including those with nuclear markers, 
we proposed the existence of only two subspecies of E. barbara (E. b. 
inserta in southern Central America, and E. b. barbara for all South 
America). All of the demographic analyses showed a very strong 
population expansion for this species in the last 400,000 YA during the 
Pleistocene. 

Keywords: Tayra, Eira barbara, mitochondrial ND5 gene, putative 
geographical subspecies, genetic diversity, phylogeography, 
population expansion during the Pleistocene 


Introduction 


Tayra {Eira barbara) is a Mustelidae (Carnivora, Mammalia) with a 
long, and slender body. Its length varies from 56 to 71 cm, not including a 
37 to 46 cm long bushy tail. Its weight ranges from 2.7 to 7 kg with males 
larger than females. This species has short, dark brown to black fur that is 
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relatively uniform across the body, limhs, and tail, except for a yellow or 
orange spot on the chest. The fur on the head and neck is much paler, 
typically tan or greyish in color. The head has small, rounded, ears, long 
whiskers and hlack eyes with a blue-green shine. The feet have toes of 
unequal length with tips that form a strongly curved line when held 
together. The claws are short and curved, but strong, being adapted for 
climbing and running rather than digging. 

This species occurs from southern Veracruz (Mexico) throughout 
Central America and across South America to northern Argentina save for 
the high Andes and the Caatinga and Cerrado (eastern Brazil; Emmons and 
Peer 1990). It is one of the most common medium-size predators 
throughout its range (Emmons and Peer 1990). 

E. barbara is a diurnal, sometimes crepuscular species (Gonzalez- 
Maya et ah, 2015), with a solitary behavior and large home range 
(Sunquist et ah, 1989). Emmons and Peer (1990) showed that the tayra 
inhabits tropical and subtropical forests, secondary rain forests, gallery 
forests, gardens, cloud forests, and dry scrub forests. Hall and Dalquest 
(1963) affirmed that it can live near human disturbed habitats. It frequently 
occurs in agricultural areas and along the edge of human settlements. Tayra 
usually inhabits areas below 1,200 m, but there are reports of it being in 
areas as high as 2,400 m (Eisenberg 1989, Emmons and Peer 1990) and it 
is common at 2,000 m (Cuaron et ah, 2016). Its diet is omnivorous, 
including fruits, carrion, small vertebrates, insects, honey and small 
vertebrates such as marsupials, rodents, and iguanids among others 
(Cabrera and Yepes 1960, Hall and Dalquest 1963, Emmons and Peer 
1990). This species is listed as Least Concern (Cuaron et ah, 2016). 

Cabrera (1957) and Hall (1981) recognized seven morphological 
subspecies, two in Central America and five in South-America: \- E. b. 
senex (Thomas in 1900). The type locality is Hacienda Tortugas, Jalapa, 
Veracruz, Mexico; 2- E. b. inserta (Allen in 1908), with the type locality in 
Ulse, Matagalpa, Nicaragua; 3- E. b. sinuensis (Humboldt in 1812), with 
the type locality for the Sinu River in the Bolivar Department in northern 
Colombia; 4- E. b. poliocephala (Traill in 1821), with type locality 
Demerara in Guyana; 5- E. b. madeirensis (Lonnberg in 1913), with type 
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locality in Humaita, Madeira River, Brazilian Amazon; 6- E. b. peruana 
(Nehring in 1886), with type locality in YuracYaku in the San Martin 
Department in Peru and 7- E. b. barbara (Linnaeus in 1758) with the type 
locality assigned hy Lonnherg (1913) to Pemanhuco, Brazil. See Figure 1. 

Although it is a relatively common species, only one preliminary study 
on its molecular population genetics and infra-specific systematics has 
been published (Ruiz-Garcia et ah, 2013). 

Therefore, we expanded upon our initial molecular population genetics 
study with mitochondrial genes of the tayra with the following main aims: 
1- To estimate the mitochondrial levels of genetic diversity in the overall 
tayra population and in some putative morphological subspecies; 2- To 
determine if there is a correlation between the molecular clades obtained in 
the phylogenetic analyses with the traditional putative morphological and 
geographical subspecies of tayras; 3- To estimate the possible temporal 
splits in the mitochondrial diversification within the evolution of the tayra; 
and 4- To determine if demographic evolutionary changes have 
characterized the natural history of the tayra. 


Materials and Methods 


Samples 

We sequenced 100 tayras at the mt ND5 gene. The samples came from 
11 countries and represent seven of the eight putative morphological 
subspecies (Table 1 & Figure 1). They are: 1- Argentina, eight individuals 
(putative E. b. barbara)', 2- Bolivia, 16 specimens (putative E. b. barbara)', 
3- Brazil, nine exemplars (four putative E. b. barbara', five putative E. b. 
madeirensis)', 4- Colombia, 12 individuals (three putative E. b. 
madeirensis', nine putative E. b. sinuensis)', 5- Ecuador, 27 specimens (four 
putative E. b. sinuensis', 23 putative E. b. madeirensis)', 6- French Guiana, 
five exemplars (putative E. b. poliocephala)', 7- Panama, one individual 
(putative E. b. inserta)', 8- Paraguay, four specimens (putative E. b. 
barbara)', 9- Peru, 17 exemplars (nine putative E. b. peruana', eight 
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putative E. b. madeirensis)] 10- Trinidad and Tobago, one individual 
(putative E. b. poliocephala). Thus, these samples represent six out of the 
seven putative morphological subspecies recognized for this species. 

The DNA of some of the tayra individuals we analyzed was extracted 
from hairs obtained from animals found alive in diverse Indian 
communities throughout Central and South America. Another fraction of 
the DNA was obtained from skins, bones, and teeth of hunted individuals 
of E. barbara. We requested permission to collect biological materials 
from these skins, bones, and teeth that were already present in the Indian 
communities. In the case of the skins, we sampled 1-2 cm^. Communities 
were visited only once. All sample donations were voluntary, and no 
financial or other incentive was offered for supplying specimens for 
analysis. For more information about sample permissions, see the 
Acknowledgment section. 



Figure 1. Map with the approximate geographical distributions of the six putative 
geographical tayra’s subspecies {Eira barbara) sequenced at the mitochondrial ND5 
gene. X represent localities where samples were obtained. 













6 Manuel Ruiz-Garcia, Nicolas Lichilm-Ortiz, Yolanda Mejia et al. 


Table 1. Samples of taya {Eira barbara), by countries, localities and 
putative geographic subspecies, sequenced at the mitochondrial 
ND5 gene for this work 


Country 

Number of 
samples studied 

Localities 

Putative geographic 
subspecies 

Panama 

1 

1 Chiriqui 

E. b. inserta 

Colombia 

12 

2 Agustrn Codazzi-Cesar 

1 PNN Los Katios-Choco 

1 Zaragoza-Antioquia 

1 Yarumal-Antioquia 

1 PNN Tama, Norte de 

Santander 

2 El Tuparro-Vichada 

1 San Martm-Meta 

1 Playa Blanca-Guainia 

1 Puerto Arara-Amazonas 

1 Leticia-Amazonas 

9 E. b. sinuensis 

3 E. b. madeirensis 

Trinidad 

& Tobago 

1 

1 Rio Claro 

E. b. poliocephala 

French 

Guiana 

5 

5 Camopi River 

E. b. poliocephala 

Ecuador 

27 

3 Yarinacocha-Pastaza 

2 Sucua-Morona Santiago 

2 La Perla-Santo Domingo 
de Tsachilas 

2 Miashi-Zamora 

2 Pillaro-Tungurahua 

2 Canelos-Pastaza 

2 Sarayaku-Pastaza 

1 Coca-Napo 

1 Misahuallr-Napo 

1 La Bonita-Napo 

1 Macuma-Morona Santiago 

1 Loreto-Napo 

1 Flushafindi-Napo 

1 Pichincha 

1 Miazal-Morona Santiago 

4 E. b. sinuensis 

22 E. b. madeirensis 
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Country 

Number of 
samples studied 

Localities 

Putative geographic 
subspecies 

Ecuador 


1 Yangana-Loja 

1 Cononaco-Pastaza 

1 Tinigua-Napo 

1 Nuevo Rocafuerte-Napo 


Peru 

17 

2 Iquitos-Loreto 

1 Caballococha-Loreto 

1 Nauta-Loreto 

1 Luceropata-Loreto 

1 Puerto Venus-Loreto 

1 Lamas-San Martin 

1 Moyobamba-San Martin 

1 Rioja-San Martin 

1 Nuevo Cajamarca-San 

Martin 

1 Bagua Grande-Amazonas 

1 Oxamarca-Cajamarca 

1 Puerto Bermudez-Pasco 

1 Bolognesi-Ucayali 

1 Seshea-Ucayali 

1 Manu-Madre de Dios 

1 Marcapata-Cusco 

8 E. b. madeirensis 

9 E. b. peruana 

Bolivia 

16 

3 Ballivian-Beni 

1 Piso Firme-Beni 

1 Nicolas Suarez-Pando 

1 Sena-Pando 

1 Franz Tamayo-La Paz 

1 Coripata-La Paz 

1 Cajuata-La Paz 

1 Totora-Cochabamba 

1 Pojo-Cochabamba 

1 Vila vila-Cochabamba 

1 Julpe River-Cochabamba 

1 El Cerro-Santa Cruz 

1 Puerto Pailas-Santa Cruz 

1 San Jose de Chuiquitos- 
Santa Cruz 

E. b. barbara 
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Table 1. (Continued) 


Country 

Number of 
samples studied 

Localities 

Putative geographic 
subspecies 

Brasil 

9 

3 Foz de Iguazu-Parana 

3 Novo Airao-Negro River- 

Amazonas 

1 Moora-Negro River- 

Amazonas 

1 Paumari-Yavari River- 

Amazonas 

1 Tres Rios-Rio de Janeiro 

4 E. b. barbara 

5 E. b. madeirensis 

Paraguay 

4 

2 Los Cedrales 

2 Loma Grande 

E. b. barbara 

Argentina 

8 

2 Salta 

1 Abra Pampa-Jujuy 

1 Humahuaca-Jujuy 

1 La Cocha-Tucuman 

1 Burruyacu-Tucuman 

1 El Dorado-Misiones 

1 San Javier-Misiones 

E. b. barbara 


Molecular Analyses 

The DNA from skins and bones was extracted using the phenol- 
chloroform procedure (Sambrook et ah, 1989), whereas DNA samples 
from hairs and teeth were extracted with 10% Chelex resin (Walsh et ah, 
1991). Primers and PCR conditions for the ND5 gene (265 bp) were 
brought to a volume of 25 pi with 13.5 pi of Mili-Q H 2 O, 3 pi of MgCb 
1 mM, 1 pi of dNTPs 0.2 mM, 1 pi of each primer (0.1 pM), 2.5 pi of 
buffer lOX, and one unity of Taq Polymerase with 50-100 ng of DNA. 
We used the primers LI 2673 and HI 2977 (5’-GGTGCAACTC 
CAAATAAAAGTA -3’ and 5^- AGAATTCTATGATGGATCATGT 3’; 
Waits et al., 1999). The PCR temperatures were 95° for 5 minutes followed 
by 10 cycles of 1 minute at 95°C, 1 minute at 64°C and 1.5 minute at 72°C, 
25 cycles of 1 minute at 95°C, 1 minute at 60°C and 1.5 minute at 72°C 
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and one final extension of 15 minutes at 72°C. All amplifications, 
including positive and negative controls, were checked in 2% agarose gels. 
Those samples that amplified were purified using membrane-binding spin 
columns (Qiagen). The PCR products were sequenced in both directions 
using the Big Dye^"^ kit in an ABI 377A automated DNA sequencer. A 
consensus of the forward and reverse sequences was determined using the 
Sequencher program. 

It is possible that some of the sequences represent numts 
(mitochondrial DNA fragments inserted into the nuclear genome) rather 
than true mtDNA (Chung and Steiper, 2008). However, we note that all 
amino acid translations of the obtained sequences showed the presence of 
initial start and terminal stop codons and the absence of premature stop 
codons. Protein translation was also checked to evaluate the possible 
presence of numts. Nevertheless, the mutations we observed were 
synonymous changes, thus suggesting that there were no numts in the 
sequences. 


Data Analyses 

Genetic Diversity 

The statistics used to determine the genetic diversity in the overall 
tayra sample and within the five South American putative tayra subspecies 
were as follows: the haplotypic diversity (Hd), the nucleotide diversity (rt), 
the average number of nucleotide differences (K), and the 0 statistic by 
sequence. These statistics were obtained using the DNAsp 5.1 software 
(Librado and Rozas, 2009). 

Phylogenetics Analyses 

The sequence alignments were carried out manually as well as with the 
DNA Alignment program (Fluxus Technology Ltd.). MrModeltest 
v2.3 software (Nylander, 2004) and Mega 6.05 software (Tamura et ah, 
2013) were applied to determine the best evolutionary mutation model. 
The Akaike and Bayesian information criteria (AIC and BIC; Akaike, 
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1974; Schwarz, 1978) were used to determine the best evolutionary 
nucleotide model in the overall sequence set of E. barbara. 

Phylogenetic trees were constructed hy using two procedures: 
Maximum Likelihood (MLT) and Bayesian analysis (BI). The ML trees 
were obtained using the RAxML v.7.2.6 software (Stamatakis, 2006). To 
select the best fitting model, 50 independent iterations were run using three 
data partitions (codon 1, codon 2, and codon 3). Additionally, 50 iterations 
were run using two data partitions (codons 1+2 combined, and codon 3). 
For each sequence data set, the GTR + G model (General Time 
Reversible + gamma distributed rate variation among sites; Tavare, 1986) 
was used to search for the ML tree and topologic support was estimated 
with 500 bootstrap replicates using GTR. 

A BI tree was completed with the BEAST v. 1.8.1 program 
(Drummond et ah, 2012). Four independent iterations were run using three 
data partitions (codon 1, codon 2, and codon 3) with six MCMC chains 
sampled every 10,000 generations for 30 million generations after a bum- 
in period of 3 million generations. We checked for convergence using 
Tracer vL6 (Rambaut et ah, 2013). We plotted the likelihood versus 
generation and estimated the effective sample size (ESS > 200) of all 
parameters across the four independent analyses to determine convergence 
and optimal results. The results from different runs were combined using 
LogCombiner vL8.0 and TreeAnnotator vL8.0 software (Rambaut and 
Drummond, 2013). A Yule speciation model and a relaxed molecular clock 
with an uncorrelated log-normal rate of distribution (Drummond et ah, 
2006) were used. Posterior probability values provide an assessment of the 
degree of support of each node on the tree. The tree was visualized in 
FigTree v. 1.4 software (Rambaut, 2012). This BI tree was used to estimate 
the time to most recent common ancestor (TMRCA) for the different 
nodes. We used a prior of 24.0 + 1 MYA (95 % confidence interval: 26.24- 
22.36 MYA) for the split between the ancestors of Eira and one 
Procyonidae, as Potos flavus. This prior followed the results of Koepfli 
et ah, (2008). 

Following Pennington and Dick (2010), the previous BI temporal 
estimates belong to one of two different approaches for inferring 
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divergence times. The first approach is based on fossil-calihrated DNA 
phytogenies. The second approach is named “borrowed molecular clocks” 
and uses direct nucleotide substitution rates inferred from other taxa. For 
this second approach, we used a median joining network (MJN) with the 
help of Network 4.6.10 software from Fluxus Technology Ltd (Bandelt 
et ah, 1999). The p statistic (Morral et ah, 1994) was estimated and 
transformed into years of divergence among the haplotypes studied. To 
determine the temporal splits, it is necessary to estimate a mutation rate at 
the mt ND5 gene. We used a nucleotide divergence of 1.22 % per each 
million years (Culver et ah, 2000), which yielded one mutation each 
309,310 years. This estimate was obtained for Felidae. In this work, we 
assumed that this mutation rate could be similar in Mustelidae. The 
networks are more appropriate for intraspecific phytogenies than tree 
algorithms because they explicitly allow for the co-existence of ancestral 
and descendant haplotypes, whereas trees treat all sequences as terminal 
taxa (Posada and Crandall, 2001). 

Heterogeneity Analyses 

Several procedures were carried out to estimate the genetic 
heterogeneity among the diverse putative tayra subspecies analyzed. To 
determine the overall genetic heterogeneity in E. barbara, we used the 
statistics Gst, yST, Nst and Fst (Nei, 1973; Hudson et ah, 1992). 
Additionally, we relied on the Hst, Kst, Kst*, Z, Z*, and Snn tests 
(Hudson, 2000), and the chi-square test on the haplotypic frequencies with 
permutation tests using 10,000 replicates to measure genetic heterogeneity. 
Also, we estimated the genetic heterogeneity by subspecies pairs within E. 
barbara. For this task, we used three procedures: 1- Exact tests with 
Markov chains, 10,000 dememorizations parameters, 20 batches, and 5,000 
iterations per batch; 2- Indirect gene flow estimates (Nm) from the Fsx 
statistic with a n-dimensional island model (Slatkin, 1985; Ruiz-Garcia, 
1993, 1994, 1997, 1999; Ruiz-Garcia and Alvarez, 2000); and 3- Kimura 
2P genetic distances (Kimura, 1980). These genetic heterogeneity statistics 
were completed with DNAsp 5.1 (Librado and Rozas, 2009) and Arlequin 
3.5.1.2 (Excoffier and Lischer, 2010). 
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Demographic Changes 

We relied on three procedures to detect possible historical population 
changes in the tayra: 1- We used the Strobeck's S statistic (Strobeck, 1987), 
Fu and Li D* and F* tests (Fu and Li, 1993), the Fu Fs statistic (Fu, 1997), 
the Tajima D test (Tajima, 1989) and the R 2 statistic (Ramos-Onsins and 
Rozas, 2002). A 95% confidence interval and probabilities were obtained 
with 10,000 coalescence permutations. 2- The mismatch distribution 
(pairwise sequence differences) was obtained following the method of 
Rogers and Harpending (1992) and Rogers et ah, (1996). We used the 
raggedness rg statistic to determine the similarity between the observed 
and the theoretical curves. 3- A Bayesian skyline plot (BSP) was obtained 
by means of the BEAST v. 1.8.1 and Tracer vL6 software. The 
Coalescent-Bayesian skyline option in the tree priors was selected with 
four steps and a piecewise-constant skyline model with 30,000,000 
generations (the first 3 million discarded as burn-in), kappa with log 
Normal [1, 1.25] and Skyline population size with uniform [0, infinite; 
initial value 80]. In the Tracer vL6, the marginal densities of temporal 
splits were analyzed and the Bayesian Skyline reconstruction option was 
selected for the trees log file. A stepwise (constant) Bayesian skyline 
variant was selected with the maximum time as the upper 95 % high 
posterior density (HPD) and the trace of the root height as the 
treeModel.rootHeight. To determine the time range for possible 
demographic changes for E. barbara, we consider that the evolution of this 
taxon occurred during the last 4 MY. 


Results 

Genetic Diversity and Phylogenetic Inferences 

The BIC showed that the best nucleotide substitution model was 
T92 -I- G (7,649.51). In contrast, the AIC detected GTR - 1 - G - 1 - 1 (5,881.29) 
as the best model. 
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The genetic diversity levels in the overall studied sample of tayra were 
very high. For the 100 individuals analyzed, we found 70 different 
haplotypes with Hd = 0.983 ± 0.006, n = 0.0422 ± 0.0048 and k = 11.175 ± 
5.117. The genetic diversity for four out of five South American putative 
morphological subspecies were very similar, all of them with very high 
genetic diversity levels (Hd = 0.991-0.960 and n - 0.0562-0.0308). The 
genetic diversity of E. b. poliocephala was somewhat lower (Hd = 0.600 
and 7t = 0.0176), although the sample size for this putative subspecies was 
the lowest (Table 2). 

The MLT and BI can be seen in Figures 2 and 3. Both phylogenetic 
trees showed that the first diverging branch represented the animal sampled 
in northcentral Panama (putatively, E. b. inserta) (MLT: low bootstrap 
28%; BI: p = 1). All of the South American specimens we analyzed were 
placed in the remaining cluster. However, although putatively animals 
from five different subspecies were included, very few significant clades 
were observed, and only partially related to the morphological subspecies. 
The first diverging cluster in the South American animals was one 
composed by three animals from northern Colombia (Cesar and Antioquia 
Departments; 50 % and p = 0.97, respectively), which corresponded with 
the putative E. b. sinuensis. Nevertheless, a Bolivian exemplar and many 
other specimens “a priori” classified as E. b. sinuensis by their 
geographical origins that did not belong to this cluster, were present in the 
BI, within this clade. Henceforth, there was only a partial correspondence 
between this clade and E. b. sinuensis. There were other interesting clades 
in both phylogenetic trees. 1- One was composed of three individuals from 
the Pacific area of trans-Andean Ecuador (61 % and p = 1, respectively), 
which also partially corresponded with E. b. sinuensis', 2- Another cluster 
was composed of individuals from different areas of Bolivia and mainly by 
individuals from northwestern Argentina (Salta, Jujuy and Tucuman 
provinces) in MLT (41 %). In the BI, this group was only composed of five 
individuals from northwestern Argentina (p = 0.96). This cluster 
was partially correlated with E. b. barbara. In the BI there was another 
cluster with several Bolivian and Argentinian specimens (p = 0.84). 
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It was separated from the first Argentinian cluster we aforementioned. 
However, as we commented for E. b. sinuensis, this relationship was 
incomplete because other individuals “a priori” classified as E. b. barbara 
were dispersed hy other clusters; 3- Another cluster of certain relevance 
was detected in the MLT and BI. It was composed of individuals of central 
Peru and one individual from the northern Peruvian Amazon (80 % and p = 
0.72, respectively). This cluster also partially supported the existence of 
the morphological subspecies E. b. peruana', 4- Small clusters of animals 
from the Ecuadorian and Colombian Amazon were present. One of them 
contained two animals from the Ecuadorian and Colombian Amazon (62 % 
and p = 1, respectively) and other three from the Ecuadorian Amazon (73 
% and p = 1; 89 % and p = 1; 28 % and p = 0.8, respectively). These very 
locally restricted clusters were inside the putative morphological 
subspecies, E. b. madereinsis. Many other individuals of E. b. madereinsis 
were distributed in clusters with other individuals “a priori” considered 
different morphological subspecies. 


Table 2. Genetic diversity in the overall sample of Eira barbara and in 
the five putative South American morphological subspecies at the mt 
ND5 gene represented by the number of haplotypes (NH), the 
haplotype diversity (Hd), the nucleotide diversity (n), and the average 
number of nucleotide differences (K) 


Eira barbara taxa 

NH 

Hd 

71 

K 

Overall Sample 

70 

0.983 +0.006 

0.0422 t0.0048 

11.175 ±5.117 

E. b. sinuensis 

12 

0.987 +0.065 

0.0562^0.0311 

14.884 ± 8.344 

E. b.poliocephala 

14 

0.600 + 0.147 

0.0176 :L0.0093 

4.667 ± 2.556 

E. b. peruana 

13 

0.989 ±0.063 

0.0517 :L0.0299 

13.703 ±7.324 

E. b. madeirensis 

26 

0.960 ±0.009 

0.0308 :l0.0071 

8.162 ±3.233 

E. b. barbara 

25 

0.991 ±0.012 

0.0412^0.0092 

10.928 ±4.111 
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Figure 2. (Continued). 
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Figure 2. Maximum likelihood tree with the 100 specimens of tayra (Eira barbara) 
sequenced at the mitochondrial ND5 gene. The number in the nodes are the bootstrap 
percentages. The procyonidae, Potos flavus, was employed as outgroup. In different 
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colors, some relevant clusters which showed a limited correspondence with some 
putative morphological geographic subspecies of E. barbara. 
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Figure 3. Bayesian tree with the 100 specimens of tayra (Eira barbara) sequenced at 
the mitochondrial ND5 gene. The three numbers in the nodes are the posteriori 
probabilities, estimated temporal splits in the nodes in millions of years, and the 95% 
high posterior density of these temporal splits. The procyonidae, Potos flavus, was 
employed as outgroup. In different colors, some relevant clusters which showed a 
limited correspondence with some putative morphological geographic subspecies of 
E. barbara. 
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Figure 4. Median Joining Network (MJN) for the haplotypes detected in 100 tayra 
{Eira barbara) sequenced at the mitochondrial ND5 gene. In light blue, one haplotype 
of E. b. inserta; in pink, haplotypes of E. b. barbara; in green, haplotypes of E. b. 
peruana; in yellow, haplotypes of E. b. madeirensis; in black, haplotypes of E. b. 
sinuensis; and in brown, haplotypes of E. b. poliocephala. Therefore, the five putative 
geographical subspecies of South American tayras and one Central America subspecies 
were represented in this analysis. Little red circles are extinct or not found haplotypes. 

The BI temporal split estimate showed an initial divergence of the 
aneestors of the Panamanian individual {E. b. inserta) and South-Ameriean 
specimens around 6.26 MYA (95%: 5.4-8.49 MYA; Miocene divergence). 
The aneestor of the elade from northern Colombia split around 3.7 MYA 
(1.55-6.37 MYA). The aneestor of the animals from northwestern 
Argentina diverged around 3.15 MYA (1.01-4.23 MYA). In contrast, the 
aneestors of animals of north-eentral Peru, trans-Andean Ecuador and one 
of the Ecuadorian Amazonas clusters diverged 2.88 MYA (1.23-4.2 
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MYA), 2.35 MYA (0.31-2.7 MYA), and 1.49 MYA (0.42-3.74 MYA), 
respectively. 

The MJN analysis revealed a view very similar to the phylogenetic 
trees (Figure 4). The major fraction of haplotypes were distributed 
irrespective of the geographical distribution of the morphological 
subspecies. For instance, the most frequent haplotypes (HI, H30, H7, H49, 
H19 and H9) included individuals of different putative subspecies: HI 
contained exemplars classified “a priori” as madereinsis, sinuensis and 
barbara', H30 and H7 were composed of madereinsis and peruana 
individuals; H49 enclosed individuals of sinuensis', H19 consisted of 
specimens of sinuensis, peruana, madereinsis, and barbara, whilst H9 
included poliocephala and peruana. Therefore, some haplotypes were 
widely distributed, which agrees quite well with extensive gene flow of 
this species across all of South America. Many of these main haplotypes 
presented other small haplotypes in star-like form, which is highly related 
to possible population expansions across the entire South American range 
of the tayra. Nevertheless, the MJN, as the phylogenetic trees, detected 
some haplotype clusters to be well delimited geographically. For example, 
there were the cases of the Central American individual (H66), the Cesar- 
Antioquia cluster (H65, H67, and H68), the trans-Andean Ecuadorian 
cluster (H45, H63, and H64), and the central Peruvian cluster (Hll, HI3, 
HI8, and H20). The MJN temporal splits were slightly less than that those 
obtained with BI, but relatively similar. The temporal splits among these 
haplotypes and the main groups can be seen in Table 3. Some of these time 
splits are interesting. The divergence between the Panamanian individual 
and H7 was estimated to occur around 4.02 + 0.31 MYA. The temporal 
divergence between clusters from the areas of northwestern Argentina, 
Cesar-Antioquia area (northern Colombia), trans-Andean Ecuador, and 
north-central Peru in reference to H7 were 3.73 + 0.65 MYA, 3.29 + 0.57, 
1.04 ± 0.31 MYA, and 2.89 ± 0.47 MYA, respectively. 

Therefore, the phylogenetics tree and the MJN analyses showed that 
the ancestor of Central and South American tayras originated during the 
Miocene or Pliocene. Also, the ancestors of some geographical groups, at 
least four of certain relevance, originated during the Pliocene and first part 
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of the Pleistocene. However, during the Pleistocene (as we will show 
later), tayra experienced a strong population expansion and many 
haplotypes expanded their geographical distributions. They superimposed 
onto the geographical areas of these older and geographical groups that 
originally differentiated during the Pliocene. 


Genetic Distances and Genetic Heterogeneity among Putative 
Morphological Subspecies of Eira barbara 

The Kimura 2P genetic distances among all of the comparison pairs of 
the six putative morphological subspecies of tayra are shown in Table 4. 
The differentiation between the Central American subspecies 
{E. b. inserta) and the five South American subspecies was elevated 
(5.5% - 7.8%), which confirmed that the Central America taxon is, at least, 
a different subspecies. It is interesting to note that the less differentiated 
South American taxon with regard to the Central American one was E. b. 
sinuensis (5.5 %). It was the South American taxon closest geographically. 
The genetic distances with the other four South American taxa ranged from 
7.I%-7.8%. 

In contrast, the genetic distances among the five South American 
subspecies were very small. They ranged from 0.1% to 1.2%. The pairs of 
subspecies with the greatest genetic distances were E. b. poliocephala-E. b. 
sinuensis (1.2%) and E. b. poliocephala-E. b. barbara (0.9%). 

The overall genetic heterogeneity for all five South American tayra 
subspecies taken together was significant (Table 5), but the genetic 
heterogeneity was relatively small. For example, the Fst and the ysx 
statistics showed values of 0.095 and 0.109, respectively. Their respective 
gene flow estimates of 4.75 and 4.10, were relatively high among the 
putative South-American subspecies. 
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Table 3. Temporal splits among different Eira barbarous lineages 
estimated by means of a Median Joining Network (MJN). 
Values of temporal splits are in millions of years. 

SD = Standard deviation 


Lineages compared 

p±SD 

Temporal 

divergence 

Between the Panamanian haplotype (inserta) 
and H49 (sinuensis) 

13.000+ 1.000 

4.021 ±0.309 

Between the Bolivian and northern-western 
Argentinian haplotypes and H7 
{madereinsis, peruana) 

12.077 + 2.097 

3.735 ±0.648 

Between the Cesar-Antioquia 
(northern Colombia) haplotypes 
and H7 {madereinsis, peruana) 

10.625 ± 1.829 

3.286 ±0.565 

Between the trans-Andean Ecuadorian haplotypes 
and H7 (madereinsis, peruana) 

3.375 + 1.008 

1.043 ±0.311 

Between the northcentral Peru and H7 
(madereinsis, peruana) 

9.333 ±1.523 

2.886 ±0.471 

Between H9 (poliocephala, peruana) and H19 
(sinuensis, peruana, madeirensis, barbara) 

1.000 + 0.500 

0.309 ±0.154 

Between H9 (poliocephala, peruana) 
and H7 (madereinsis, peruana) 

1.363 +0.454 

0.422 ±0.140 

Between H19 (sinuensis, peruana, madeirensis, 
barbara) and H7 (madereinsis, peruana) 

0.454 + 0.454 

0.141 ±0.141 

Between H49 (simuensis) and H7 
(madereinsis, peruana) 

1.429 + 0.714 

0.441 ±0.220 

Between H30 (sinuensis) and H7 
(madereinsis, peruana) 

0.625 ± 0.625 

0.193 ±0.193 

Between H9 (poliocephala, peruana) 
and HI (madeirensis, sinuensis) 

2.400 + 0.600 

0.742 ±0.185 

Between H19 (sinuensis, peruana, madeirensis, 
barbara) and HI (madeirensis, sinuensis) 

1.200 ±0.600 

0.371 ±0.185 

Between H49 (sinuensis) and HI (madeirensis, 
sinuensis) 

2.454 ±0.818 

0.759 ±0.253 

Between H7 (madereinsis, peruana) 
and HI (madeirensis, sinuensis) 

0.643 ± 0.643 

0.198 ±0.198 

Between H30 (sinuensis) and HI 
(madeirensis, sinuensis) 

1.500 ±0.790 

0.463 ±0.231 
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Table 4. Kimura 2P genetic distance (Kimura, 1980) in percentages 
(%) among six different morphological subspecies of Eira barbara 
(Mustelidae) (below main diagonal) and standard deviations in 
percentages (%) (above main diagonal) at the mt ND5 gene. 

\ = E.b. barbara', 2 = E. b. peruana', 3 = E. b. madeirensis', 

4 = E. b. poliocephala', 5 = E. b. sinuensis', 6 = E. b. inserta', 

7 = Potos flavus (Procyonidae) 



1 

2 

3 

4 

5 

6 

7 

1 


0.1 

0.1 

0.4 

0.1 

1.4 

3.5 

2 

0.2 


0.1 

0.2 

0.1 

1.5 

3.5 

3 

0.2 

0.1 


0.4 

0.1 

1.5 

3.7 

4 

0.9 

0.5 

0.8 


0.5 

1.5 

3.2 

5 

0.3 

0.5 

0.4 

1.2 


1.2 

3.5 

6 

7.1 

7.7 

7.4 

7.8 

5.5 


3.5 

7 

27.5 

27.5 

28.5 

24.6 

26.5 

24.9 



Table 5. Overall genetic heterogeneity and gene flow (Nm) statistics 
for the flve putative South American subspecies of Eira barbara at 
the mt ND5 gene. * P < 0.05; ** P < 0.01 


Estimated Genetic 

Differentiation 

P 

Gene flow 

313.351 df=272 

0.0429* 

Gst = 0.0336 

Nm = 14.40 

Hst= 0.0236 

0.0001** 

ysT = 0.0951 

Nm = 4.75 

Kst = 0.0559 

0.0001** 

Nst= 0.1108 

Nm = 4.01 

Kst* = 0.0362 

0.0001** 

Fst = 0.1086 

Nm = 4.10 

Zs = 2279.890 

0.0001** 



Zs* = 7.360 

0.0001** 


S„„ = 0.511 

0.0001** 



The analysis of subspecies pair comparisons with exact prohahility 
tests (Table 6) only showed two significant pairs: E. b. madereinsis-E. b. 
barbara (p = 0.0066 + 0.0017) and E. b. madereinsis-E. b. poliocephala 
(p = 0.0135 + 0.0034). In this analysis, the Central American taxon was 
deleted because only one sequence was analyzed. The estimates of Nm by 
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subspecies pair comparisons (Table 7) clearly yielded that the values lower 
than 1 (which is considered the limit for low gene flow; Wright, 1943) 
always implied the Central American taxon: E. b. inserta-E. b. peruana 
(Nm = 0.584), E. b. inserta-E. b. barbara (Nm = 0.459), E. b. inserta-E. b. 
madereinsis (Nm = 0.286), E. b. inserta-E. b. poliocephala (Nm = 0.137). 
Only a South American taxon, E. b. sinuensis, (Nm = 1.202) had a more 
substantial gene flow with E. b. inserta. This agrees quite well with that 
determined in the phylogenetic trees and in the genetic distance analysis. 
The gene flow estimates among the South American taxa were all above 1, 
ranging from 2.469 to 11.847. These values strongly correlate to elevated 
historical gene flows among the populations of tayra throughout South 
America. 


Table 6. Exact probability tests (P) (below main diagonal) and 
standard deviations (above main diagonal) among six different 
putative morphological subspecies of Eira barbara by means of 
the mt ND5 gene. \ = E.b. barbara', 2 = E. b. peruana', 

3 = E. b. madeirensis', 4 = E. b. poliocephala', 5 = E. b. sinuensis', 
6 = E. b. inserta. * = Significant probability 



1 

2 

3 

4 

5 

6 

1 


0.0000 

0.0017 

0.0138 

0.0197 


2 

1.0000 


0.0071 

0.0132 

0.0088 


3 

0 . 0066 * 

0.0679 


0.0034 

0.0201 


4 

0.2307 

0.3472 

0 . 0135 * 


0.0036 


5 

0.4044 

0.4796 

0.1032 

0.0893 



6 













Demographic Evolutionary Changes in the Tayra 

All of the demographic change statistics indicated population 
expansion in the tayra (Strobeck's S statistic, P = 0.0001; Tajima D = - 
2.258, p = 0.0040; Fu & Li D* = -3.175, p = 0.0115; Fu & Li F* = -3.335, 
P = 0.0052; Fu’s Fs = -52.632, P = 0.00001; and Rz = 0.037, P = 0.0041). 
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Table 7. Gene flow (Nm) estimates (below main diagonal) among six 
different putative morphological subspecies of Eira barbara by means 
of the mt ND5 gene. 1 = E. b. barbara', 2 = E. b. peruana', 3 = E. b. 
madeirensis', 4 = E. b. poliocephala', 5 = E. b. sinuensis', 6 = E. b. inserta 



1 

2 

3 

4 

5 

6 

1 







2 

9.8245 






3 

9.2893 

11.8471 





4 

2.6879 

5.3236 

2.2686 




5 

7.1919 

5.8729 

4.8435 

2.4691 



6 

0.4597 

0.5843 

0.2862 

0.1373 

1.2019 




Figure 5. Historical demographic analysis by means of the mismatch distribution 
procedure (pairwise sequence differences) for the mitochondrial ND5 gene studied in 
the overall sample of Eira barbara. The analysis showed a clear population expansion 
of this species during the Pleistocene. 

The mismatch distributions also indicated population expansion (rg - 
0.0040, P = 0.00280) (Figure 5). Assuming one year as one generation in 
the tayra, the population expansion began 343,586 YA, during the 
Pleistocene. 
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Time 


Figure 6. Bayesian skyline plot analysis (BSP) to determine possible demographic 
changes across the natural history of the overall sample of tayra {E. barbara) 
sequenced at the mitochondrial ND5 gene. The analysis showed a clear female 
population expansion during the Pleistocene. On the x-axis, time in millions of years; 
on the y-axis, log effective population size of females. 

The BSP analyses also determined a strong female population 
expansion during the Pleistocene for the tayra (Figure 6). The analyses 
showed the beginning of the expansion around 400,000 YA, very similar to 
the temporal estimate previously showed. Therefore, there is 
incontrovertihle evidence that the tayra experienced a strong population 
expansion during the Pleistocene, as was previously suggested hy the MJN 
analysis. 


Discussion 


Genetic Diversity 

The levels of nucleotide diversity found in E. barbara (n - 0.0422) 
were quite high. They were higher than those found in many other 
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neotropical carnivores. For example, they were higher than three fox 
species (Lycalopex culpaeus: n - 0.008, Ruiz-Garcra et ah, 2013a; 
Lycalopex sechurae: n — 0.015, Ruiz-Garcra et ah, 2013a; Cerdocyon 
thous: n - 0.019, Tchaika et ah, 2006), three otter species (Lontra felina: 
n - 0.005, Vianna et ah, 2010; Lontra longicaudis: n - 0.011, Ruiz-Garcra 
et ah, 2017a; Pteronura brasiliensis: n - 0.0086, Ruiz-Garcra et ah, 
2017a), and two vulnerable Neotropical cats (Leopardus jacobita: n - 
0.0047, Ruiz-Garcra et ah, 2013b; Leopardus guigna: n - 0.00461, 
Napolitano et ah, 2013). The values of E. barbara were similar to those 
found in certain Neotropical cats, which are characterized by very elevated 
genetic diversity levels {Puma yaguaroundi: n - 0.0661; Ruiz-Garcra and 
Pinedo-Castro, 2013 and Ruiz-Garcra et ah, 2017b; Leopardus pajeros: n - 
0.0513, Ruiz-Garcra et ah, 2013b; Leopardus pardalis: ti - 0.068, Eizirik 
et ah, 1998; Leopardus wiedii: n - 0.035-0.074, Ruiz-Garcra et ah, 2017c). 

These high levels of genetic diversity in “a priori” neutral markers, as 
that we studied, could be related with the fact that many other genes in the 
genome of a given species contains enough variability for the action of 
natural selection (Kimura, 1986). This could be the origin of the great 
morphological variation and behavior plasticity found in the tayra. 
Emmons and Ereer (1990) determined that tayras could live in a wide 
variety of habitats such as tropical and subtropical forests, primary rain 
forest (as throughout the Amazonian forest in Brazil, Peru, Colombia and 
Ecuador), secondary rain and gallery forests (as in the Llanos of Venezuela 
and Colombia), gardens, plantations, and cloud forests. They also inhabit 
dry scrub and deciduous forests (as in the Pantanal in Brazil, Paraguay and 
Bolivia) and tall grass savannas (as in Argentina, Bolivia, and Paraguay). 
Sunquist et ah, (1989) showed that the extreme plasticity of this species for 
habitat preferences, activity periods, and diet preferences may reduce 
interspecific competition between E. barbara and other carnivores. This 
could also be the explanation why Konecny (1989) found no significant 
habitat preference for this mustelid in Belize. The abundance of the tayra 
throughout much of Central and South America may be a consequence of 
its ecological flexibility compared to sympatric carnivores. Associated 



30 Manuel Ruiz-Garcia, Nicolas Lichilm-Ortiz, Yolanda Mejia et al. 


with this, the tayra is a generalist predator, consuming a variety of fruits, 
carrion, small and medium vertebrates, insects, and honey (Cahrera and 
Yepes, 1960; Galef et al., 1976; Hall and Dalquest, 1963; Konecny, 1989; 
Sunquist et ah, 1989). 


Genetic Heterogeneity, Gene Flow and the Systematics of 
the Tayra 

Our results clearly showed that the specimen sampled in Central 
America was highly divergent from all of the individuals sampled in South 
America. However, the genetic heterogeneity among the putative 
morphological subspecies of South American tayras, although significant, 
is very small as we found with the Fst statistic, exact probability tests, and 
genetic distances. The indirect gene flow estimates were clearly higher 
than 1. Wright (1943) stated that in an island model, if Nm > 1, then gene 
flow is important enough to erase the genetic heterogeneity among 
populations. In a stepping-stone model, this amount must be larger than 4 
(Trexler, 1988). In both models, Nm < 0.5 means that the populations are 
highly disconnected from a reproductive point of view. For instance, the 
gene flow estimates between E. b. inserta and E. b. poliocephala 
(Nm = 0.137), E. b. inserta and E. b. madereinsis (Nm = 0.286) and E. b. 
inserta and E. b. peruana (Nm = 0.459) showed that the Central American 
taxon is completely isolated from these three South American taxa. 
However, recall that certain genetic relatedness was detected between the 
Central American taxon and the most northern South American taxon (E. 
b. sinuensis). Additionally, we found several gene flow comparison pairs 
between South American taxa, such as E. b. madereinsis-E. b. barbara 
(Nm = 9.289) and E. b. madereinsis-E. b. poliocephala (Nm = 2.168). 
They were elevated, although these comparison pairs showed significant 
heterogeneity. This might be explained according to Alledorf and Phelps 
(1981), who argued that the most correct interpretation of Nm > 1 is that 
the populations share the same alleles, although not necessarily with the 
same allele frequencies. By means of simulations, these authors showed 
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that significant allele divergence occurred in 50% of the generations with a 
gene flow of Nm = 50. Significant allele divergence happened on most 
occasions when Nm =10. 

The tayra seems to have a strong dispersion capacity, which could he 
related with these high gene flow estimates detected for all the putative 
South American subspecies. For instance, in the Venezuelan and 
Colombian Llanos, tayras are usually found along gallery forests. 
However, tayras cross these extensive grasslands at night, presumably 
moving from one forest to another covering long distances (Defier, 1980). 

Taking into consideration all these facts, we suggest that the six 
putative morphological subspecies analyzed could be reduced to two 
different subspecies: E. b. inserta for southern Central America and E. b. 
barbara for all South America. The name should be E. b. barbara because 
it was given by Linneus in 1758 versus E. b. sinuensis given by Humboldt 
in 1812, E. b. peruana given by Nehring in 1886, E. b. madeirenis given 
by Lonnberg in 1913 and E. b. poliocephala given by Traill in 1821. In 
reference to the putative northern Central America subspecies (E. b. 
senex), we cannot add any comment on its systematics because no 
individual of this putative subspecies was analyzed. Therefore, it is 
essential to sample tayras from this putative subspecies to determine its 
relationships with other tayra taxa. 

Here we suggest another alternative point of view. As we showed, the 
first tayra’s splits originated during the Miocene-Pliocene and beginning of 
the Pleistocene. However, during the Pleistocene, tayra experienced a 
strong population expansion and many haplotypes expanded their 
geographical distributions and they became superimposed on the 
geographical areas of older geographical groups that originally 
differentiated during the Pliocene. We suggest that future studies analyze 
nuclear genes to determine if there was hybridization between the older 
geographical groups (northern Colombia, part of Bolivia and northwestern 
Argentina, trans-Andean area of Ecuador, and north-central Peru) and the 
tayra’s population that expanded throughout South America during the 
Pleistocene. If data support this, then our view of a unique tayra’s 
subspecies in South America should be valid. On the contrary, if there was 
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little or no hybridization between the original groups and Pleistocene 
colonizers in sympatry, then the number of subspecies in South America 
could be higher. Therefore, the northern Colombian population (Cesar, 
Antioquia and possibly nearby areas) should be named E. b. sinuensis and 
the northern central Peruvian population should be named E. b. peruana. 
Also, the Bolivian and especially the northwestern Argentinian population 
should be defined as a new subspecies (tentatively E. b. saltensis). The 
trans-Andean Ecuadorian population should be defined as a new 
subspecies (tentatively E. b. aequatorialis). The remaining populations of 
tayra in South America should be named as E. b. barbara. Additionally, 
the range distribution of E. b. sinuensis and E. b. peruana should be more 
restricted than traditionally accepted (see Presley, 2000). Even, if some 
reproductive isolation mechanism had emerged between the Central and 
the South America tayras due to the old split estimated during the 
Miocene-Pliocene, both tayra populations should be consider two different 
species (E. barbara and E. inserta', this last should be E. senex if both 
Central American forms of tayra were genetically undifferentiated because 
senex was firstly named by Thomas in 1900 and inserta was named by 
Allen in 1908). 

Only future nuclear genetic studies can clarify which of the two points 
of view is more acceptable. 


Temporal Splits in the Tayras 

Our results showed that the divergence between the Central and the 
South American tayras occurred around 6.3 to 4 MYA (Miocene or 
Pliocene periods depending of the temporal estimation). Johnson and 
O'Brien (1997) and Johnson et ah, (2006) showed that seven of the eight 
primary lineages of felids radiated in the early part of the Late Miocene 
(10.8-6.2 MYA). There was a noteworthy cooling of the global climate 
near the end of the Middle Miocene. This period of cooling coincides with 
formation of a permanent Antarctic ice sheet in the Middle and in the Late 
Miocene and an Arctic ice sheet in the Pliocene. A large peak of 
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diversification in many vertebrate taxa occurred during the Pliocene epoch. 
The cold and dry climate during the Pliocene, coincides with the onset of 
high latitude glacial cycles, causing an explosive expansion of low- 
biomass vegetation, including grasslands and steppe at mid-latitudes and 
development of taiga at high latitudes of Eurasia and North America. 
These changes were correlated with the diversification of prey species such 
as muroid rodents and passerine birds that exploited these new habitats, 
which in turn provided new niches for little or medium carnivores, such as 
the tayra. Additionally, this Pliocene period agrees quite well with the last 
phase of rising of the Andes as shown by Dollfus (1974) and Clapperton 
(1993) (see, for instance, the rising of the “tablazos” of Piura, Peru) and 
very high volcanic activity in the Andes, with replacement of rainforests 
with steppe and grassland environments. 

Our initial divergence estimates in the tayra agrees relatively well with 
other molecular studies in the reconstruction of the phylogenetic 
relationships in the Mustelidae (Bryant et ah, 1993; Dragoo and Honeycutt, 
1997; Koepfli and Wayne, 1998, 2003; Sato et ah, 2003, 2004; Flynn et ah, 
2005; Fulton and Strobeck, 2006; Koepfli et ah, 2008). Two molecular 
studies are fundamental to understanding the phylogenetics of the 
Mustelidae (Koepfli and Wayne, 2003; Koepfli et ah, 2008). In the first 
study, the authors used five nuclear gene segments and the mt Cyt-b gene. 
The genes APOB, FES, GHR, RHOl and mtCyt-b clustered E. barbara 
together with Martes americana, Martes pennanti and Gulo gulo. On the 
other hand, CHRNAl clustered E. barbara with Meles meles and Arctonyx 
collaris. The major part of the trees generated by these authors showed that 
Mustelinae and Melinae were polyphyletic within the Mustelidae, whereas 
Lutrinae was monophyletic. The authors of the second study analyzed 22 
nuclear and mitochondrial gene segments and determined Mustelidae to 
consist of seven primary groups. These groups include four major clades 
and three monotypic lineages. It also included Eira barbara clustered into 
the subfamily Martinae, together with Martes and Gulo, the most divergent 
taxa within this subfamily (100 % of bootstrap and posterior probability of 
1). In that study, the branch of E. barbara diverged from the other 
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Mustelidae around 6.7-7.7 MYA (calculated using the mean values), which 
agrees with our estimate. 

These molecular results are not in conflict with the fossil record we 
know and understand for Mustelidae in America. Many species of 
Mustelidae appeared in North America during the Late Miocene-Early 
Pliocene. For instance, Cemictis hesperus, from the Pinole Tuff Local 
Fauna of California, has been dated radiometrically to have lived 5.3-5.5 
MYA (Tedford et ah, 2004; Baskin, 2011). Other cases of extinct genera 
appearances are Trogonictis and Sminthosinis as well as extant genera such 
as Lutra and Mustela during the Hemphilian period (4.7-5.9 MYA, 
Tedford et ah, 2004). There is also Legionarictis fortidens from the 
Barstovian (Middle Miocene) marine Temblor Formation in California 
(Tseng et ah, 2009). This form has shown a very close resemblance to 
other Miocene Mustelidae genera such as Dehmictis, Eirictis, Iberictis, and 
Trochictis, all from the Old World. The form also closely resembles 
Sminthosinis, Trigonictis (all from the New World), and, especially, with 
two extant genera, Galictis and Eira (Ray et ah, 1981; Ginsburg and 
Morales, 1992; Baskin, 1998; Ginsburg, 1999). In fact, the cladistic 
analysis of Tseng et ah, (2009) determined that this genus could be an 
evolutionary basal stage closely related to Eira. At the Longdan Fauna of 
the Gansu Province in China, a similar fossil to Eira was found by Qiu 
et ah, (2004) from the Late Pliocene {Eirictis robusta', 2.58-2.15 MYA). 
However, the cladistic analysis of Tseng et ah, (2009) determined that this 
mustelid was not very close to Eira and it is probably younger than Eira. 
Two possible fossil species of Eira have been described from post- 
Pliocene deposits of Maryland and Virginia under the names Galera 
macrodon and G. perdicida. However, the former has been assigned to 
Trigonictis based on additional material collected from deposits of the 
Blancan land mammal age from Washington, Idaho, Nebraska, Kansas, 
Texas, North Carolina, and Florida (Ray et ah, 1981). Trigonictis is 
considered an intermediate form between Galictis and Eira and could be 
ancestral to both. The second species may be Mephitis (Alston, 1882). In 
addition, extinct species of Eira were noted from the Pliocene of the 
Eastern Hemisphere, but specific names were not given (Scott, 1937). 
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Hershkovitz, (1972) claimed that Eira and other endemic monotypic 
mustelid genera, such as Lyncodon and Pteronura, may have evolved in 
South America and moved north as part of the north and south American 
interchange across the Panamanian land bridge. In contrast, fossil records 
suggest that Eira may have a North American origin (Ray et ah, 1981). 

Therefore, the process of diversification within Eira could he at the 
end of the first mustelid diversification peak or at the beginning of the 
second mustelid diversification peak detected by Koepfli et ah, (2008). In 
either case, this mitochondrial diversification process occurred before the 
Panamanian land bridge (3-2.5 MYA). Therefore, Eira's could have 
radiated in North America before South America in concordance with the 
fossil record (Ray et ah, 1981) and against the view of Hershkovitz (1972). 
If so, tayras arrived in South America before the complete formation of the 
Panamanian land bridge, coinciding with the Choco-Panama island bridge 
(Galvis, 1980), which could have been used by the ancestors of the current 
E. barbara to colonize northern South America from Central America. 
During the upper Pliocene orogeny, the present Tuira, Atrato and Sinu 
river basins and the nearby lowlands were raised above sea level. Thus, the 
mountains of southern Central America and of the northern Andes were 
uplifted to about their present elevation (Van der Hammen, 1961). 
Although the Nicaraguan, Panamanian and Colombian portals remained 
open (upper Miocene-middle Pliocene), numerous volcanic islands existed 
from the lower Atrato Valley and the Tuira River Basin of eastern Panama 
to the Nicaraguan portal. They could have been used by Eira's ancestor to 
migrate southward. The Cuchillo Bridge of the Uraba region, connecting 
the Tertiary Western Colombian Andes with the Panamanian islands was 
probably above sea level during this period. Henceforth, tayras could be 
another “island hopper” species (Simpson 1950, 1965, 1980). 

Our results could also be considered as indirect evidence of a Miocene 
origin of the Isthmus of Panama (Montes et ah, 2012, 2015). Indeed, the 
Isthmus of Panama formation began earlier and seems to be associated 
with the Northern Andean uplift, around 24 MYA (Farris et ah, 2011). 

Therefore, the tayra could have arrived in South America before other 
Mustelidae. Koepfli et ah, (2008) claimed that genera and species of 
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mustelids found in South America today are largely descended from North 
American immigrants that arrived as part of the GABI following the rise of 
the Panamanian isthmus, 3.0-2.5 MYA. Several informational points 
support the statement of Koepfli et ah, (2008). 1-For example, there is the 
clade of New World otters where Lutra canadensis is sister to Lontra 
felina, and Lontra longicaudis. The latter two species are found in South 
America and are estimated to have diverged 2.8-3.4 MYA (95% HPD: 
1.6-5.2 MYA) overlapping with the formation of the Panamanian land 
bridge. 2-The long-tailed weasel, Mustela frenata, ranges from North 
America to northern South America. In addition, Mustela africana and 
Mustela felipei are endemic to South America. Fossil evidence clearly 
indicates that Mustela colonized South America from the north, apparently 
well after the Panamanian isthmus was in place. 3- Fossils of the current 
Lyncodon patagonicus, (and a fossil form very related as Lycodon bosei) 
and Stipanicicia sp (closely related to the extant Galictis cuja) have been 
registered in the Ensenadean and Bonaerense periods of the Argentinian 
Pleistocene (Forasiepi et ah, 2007). However our results could be ratified 
by other results provided by the same authors (Koepfli et ah, 2008). They 
found that Pteronura, Galictis and Lira could have a Eurasian origin for 
each genus with posterior diversification in North America. Eor example, 
Pteronura may be related to the extinct genus Satherium from the Pliocene 
of North America. Additionally, Lira may be related to Trigonictis and 
Legionarictis, also North American fossils (Tseng et ah, 2009). Eossil 
evidence suggests mustelids colonized the New World across Beringia 
during different intervals when the land bridge between Eurasia and North 
America was open. Multiple genera of mustelids migrated into North 
America during the Late Miocene (around 11.2-5.3 MYA), prior to the first 
opening of the Bering Strait 5.4-5.5 MYA, which severed the route across 
Beringia. Many genera that colonized North America during the Late 
Miocene or earliest Pliocene became extinct. However, Lira could be one 
of the surviving genera that began to diversify in North America and also 
in South America if we accept that they arrived in South America before 
the complete formation of the Panamanian land bridge. 
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The mitochondrial diversification of the oldest groups of South 
American tayra occurred 3.7-2.3 MYA. This coincided with the climatic 
changes that originated from the completion of the Panamanian land bridge 
(3.1-2.8 MYA; Marshall et ah, 1979, 1982; Marshall, 1985, 1988; Wehh, 
1985, 1997; Coates and Ohando, 1996) in the Last Pliocene. 
Diversification in the tayra occurred close to the Gelasian period (2.5-1.8 
MYA), a period characterized hy the last stages of a global cooling trend 
that led to the quaternary ice ages (International Commission on 
Stratigraphy 2007). Around 2.5 MYA, the Andean forests were 
transformed into open cold dry savannah (‘paramo’), which could have 
potentially isolated populations of different species. They could have 
crossed the Northern Andes coming from Central America. Van der 
Hammen (1992) demonstrated that the mean temperature in the Colombian 
Andes was 4 ° C lower than today. He also stated that the rain level 
descended below the level reported for today (500-1,000 mm). At 2,500 
meters above sea level (masl), the temperature was 10 °C lower than it is 
today. Tayra’s diversification could have been affected by the rapid uplift 
that resulted in a significant elevation of the Northern Andes. The 
mountain range’s height climaxed around 2.7 MYA when the northern 
Colombian Andes reached its present day elevation (Gregory-Wodzicki, 
2000). This also coincides with the last formation of the Central Andes. All 
of the Andean chain between Cajamarca and Huancavelica in Peru 
appeared by volcanism in this period. 

Much of the mitochondrial diversification process of the typical 
Pleistocene colonizer haplotypes occurred around 1.5-0.8 MYA. This 
divergence could have been initiated by the pre-Pastonian glacial period 
(1.3-0.8 MYA), which had the highest glacial peak of the first Quaternary 
glacial period (Giinz). This glacial period was extremely dry, and there was 
a great degree of forest fragmentation. This period was a time for 
haplotype diversification. It was also a time of separation for many 
carnivores as it was previously determined for the Pampas cat (Cossios 
et al. 2009), and for the foxes of the Lycalopex genus (Ruiz-Garcia et al. 
2013a). Around 1.3 MYA, the Buenos Aires’s fauna transformed into a 
typical semi-arid Patagonian fauna, represented by the guanaco. 



38 Manuel Ruiz-Garcia, Nicolas Lichilm-Ortiz, Yolanda Mejia et al. 


Lestodelphys and Lyncodon. Therefore, the climate was considerably 
colder and drier than today and could have influenced the mitochondrial 
fragmentation within the tayra. 

The strong population expansion detected around 0.4 MYA for the 
tayra agrees well with an interglacial period (0.39-0.20 MYA, West 1967) 
characterized hy higher temperatures and humidity and forest expansions 
(Hoxniense in the British Islands, Yarmouth in North America, Holstein in 
northern Europe and Mindel-Riss interglacial period in central Europe). 

Euture analyses with nuclear markers are needed as well as samples 
from Central America (especially, southern Mexico, Guatemala, Belize 
and Honduras), the Pacific areas of Colombia and Ecuador, other areas of 
Central Peru, and the Guyana shield. These markers and additional samples 
will help us to determine the exact number of subspecies or ESUs (Moritz, 
1994). This information is crucial for the development of effective 
conservation plans for this species. 
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